•  • 


security  classification  of  this  page  f**s«  Qat*  twtfd) 


|  REPORT  DOCUMENTATION  PAGE 

READ  INSTRUCTIONS 

BEFORE  COMPLETING  rORM 

1.  RCFORT  nunecn 

Technical  Note  BN-1131 

2.  GOVT  ACCESSION  NO^ 

1.  RECIPIENT'S  CATALOG  HUMBER 

4.  TITLE  fan*  Svktnlt) 

BASIC  PROBLEMS  OF  A-POSTERIORI  ERROR  ESTIMATION 

Filial  Life  of  Contract 

s.  performing  org.  report  number 

7.  AU  TMORf  •> 

1  2  3 

I.  Babuska  -  L.  Plank  -  R.  Rodriguez 

s.  contract  or  grant  Humaino 

*0NR  N00014-90-J-1090 
'CCR -88-20279  (NSF) 

^CCR -88-20279  (NSF) 

s.  performing  organization  name  and  AOORCSS 

Institute  for  Physical  Science  and  Technology 
University  of  Maryland 

College  Park,  MD  20742-2431 

10.  PROORAM  ELEMENT.  PROJECT.  TASK 
AREA  0  WORK  UNIT  NUMBERS 

II.  CONTROLLING  OFFICE  NAME  ANO  AOORESS 

Department  of  the  Navy 

Office  of  Naval  Research 

12.  REPORT  OATS 

February  1992 

12.  NUMBER  OF  PAGES 

28 

"UTMONltORING  AGENCY  NAME  A  AODRESSfll  GUtml  AM  CentrellAiJ  01  tin) 

IS.  SECURITY  CLASS.  (•!  reperlj 

It  a.  DECLASSIFICATION/  DOWNGRADING 
SCHEDULE 

is.  distribution  statement  f#i  mi*  s««u 

Approved  for  public  release:  distribution  unlimited 

17.  DISTRIBUTION  STATEMENT  fal  IA*  atctracl  nlm4  Im  *!•«*  7*.  If  *11  •ml  Am  *M») 

IS.  SUPPLEMENTARY  NOTES 

IS.  KEY  WORDS  fCatiHmie  «n  nwmn  «l*  II  hkhint  fM 

Abstract.  This  paper  addresses  b 

a-posterlorl  error  estimation.  I 

Until tr  *F  S fee*  mein) 

aslc  questions  regard: 

t  analyzes  In  detail  t 

triangular  meshes.  It 

n  the  geometry  of  the 

lutlon.  The  theoretic 

ur  nuaerlcal  experlaer 

ng  a  one  parametric  fa 

daptlvely  constructed 

ng  reliability  of 

in  a-posterlorl 

ZO.  ABSTRACT  <  estimator  /or  linear  elements  on 

the  effectlvlty  Index  depending  o 

smoothness  of  the  approxlaated  so 

concordance  with  the  results  of  o 

analyzes  In  one  dlaenslonal  settl 

and  their  optlaal  selection  for  a 

gives  the  bounds  of 

triangles  and 

:al  results  are  In 

its.  The  second  part 

tally  of  the  estlaators 

meshes. 

DO  ,  j  ah*?)  1473  COITION  OF  I  NOV  ••  IS  obsolete 

S/N  010I- LMP4»  4401  SECURITY  CLASSIFICATION  OF  TniS  RAO*  (■*•>  Out*  InVared) 


a.-  eor 

'  uXAftl 
(  ’  Ti.6 

:ntfvui.."«i«»a  j 

.7  isti  reflation . 

1  ■/ - 

In/ 

A.’ftUa»liity  ;0i«6 
Uva'l  aad/ar 
Dlst  i  Spoalal 

IM 


BASIC  PROBLEMS  OF  A-POSTERIORI  ERROR  ESTIMATION 


I.  Babuska , ^  L.  Plank , ^  R.  Rodriguez"* 


r*v' 


Abstract.  This  paper  addresses  basic  questions  regarding  reliability  of 
a-posteriori  error  estimation.  It  analyzes  in  detail  an  a-posteriori 
estimator  for  linear  elements  on  triangular  meshes.  It  gives  the  bounds  of 
the  effectivity  index  depending  on  the  geometry  of  the  triangles  and 
smoothness  of  the  approximated  solution.  The  theoretical  results  are  in 
concordance  with  the  results  of  our  numerical  experiments.  The  second  part 
analyzes  in  one  dimensional  setting  a  one  parametric  family  of  the  estimators 
and  their  optimal  selection  for  adaptively  constructed  meshes. 
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1.  INTRODUCTION 


During  the  last  ten  years,  starting  essentially  with  [2-5],  great  progress 
has  been  made  in  the  theory  and  practice  of  a-posteriori  error  estimation. 
A-posteriori  error  estimation  is  then  directly  related  to  adaptive 
construction  of  the  meshes.  For  survey  of  various  results,  we  refer  to  [6], 
[7],  [15],  [24]  and  the  references  therein.  Many  estimators  were  developed 
and  analyzed  either  by  numerical  experimentation,  mathematical  theory,  or 
both.  Nevertheless,  there  is  no  sufficiently  detailed  mathematical 
understanding  of  various  estimators,  especially  regarding  their  reliability 
(e.g.,  bounds  of  the  effectivity  indices,  etc.). 

In  this  paper,  we  analyze  in  detail  a  particular  a-posteriori  error 
estimator  for  linear  triangular  elements.  We  show  the  dependence  of  its 
effectivity  index  on  various  factors  such  as  geometry  and  smoothness  of  the 
solution.  We  present  theoretical  results  and  numerical  benchmark 
computations.  Further,  we  analyze  a  one  parametric  family  of  a-posteriori 
error  estimators  in  one  dimensional  setting.  We  show  that  for  adaptively 
constructed  meshes  (which  are  equilibrated),  asymptotic  exactness  of  the 
estimator  can  be  achieved  for  unsmooth  solutions  by  an  optimal  choice  of  the 
parameter. 


2.  BASIC  NOTIONS  IN  THE  THEORY  OF  THE  A-POSTERIORI  ERROR  ESTIMATION  FOR  A 
MODEL  PROBLEM 

Let  £1  c  IR2  be  a  polygon  with  the  boundary  T.  Denote  by  H  =  Hq  (ST2 ) 
the  standard  Sobolev  space,  and  consider  the  bilinear  form  B(u,v)  on  H  x  H: 


with 


B(u, v) 


f  f£u  3v 

Jn  [ax  ax 


+ 


3u  av] 

ay  3yJ 


dxdy. 


[|uj|E  =  (B(u,  v)  )1/2. 


Let  f  €  L2(0)  and 
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F(v)  =  f  fv  dxdy  €  H' 

Jn 

be  a  functional  on  H.  Then,  by  uq  e  H  we  will  denote  the  (exact)  solution 
of  the  problem  T 

(2.1)  B ( uq , v )  =  F(v) ,  V  v  e  H. 

Obviously,  uq  is  uniquely  defined  by  (2.1) 

Let  M  be  a  class  of  triangularizations  of  £1;  by  7  €  M  we  denote  the 
triangular  mesh  and  by  t  -  7  the  triangular  element.  Let  S(7)  c  H  be 
the  set  of  all  functions  which  are  linear  on  every  element  t  €  7.  Then,  as 
usual,  by  ug  €  S(7)  we  denote  the  finite  element  solution  of  the  problem  T: 

(2.1)  B(u  , v )  =  B(u  , v )  =  F(v) ,  V  v  e  S(J). 

S  0 

We  associate  to  every  t  €  7  an  error  indicator  t)(t)  £  0,  which  value 
depends  only  on  ug  and  f  on  t  and  on  the  neighboring  elements.  The 
concrete  form  of  some  error  indicators  will  be  given  below. 

We  define  the  error  estimator 

(2.2)  g(7)  =  ^  £  tj2(t) 

xe7 

We  will  be  interested  in  estimating  the  error  in  the  energy  norm,  i.e.  we 
assume  that  6(7)  =  flu  -u  II  .  Further  we  define  the  effectivity  index  0(7) 

11  0  S  1 E 

of  the  estimator  8: 

(2.3)  0(7)  = 

The  error  Indicator,  estimator  and  effectivity  index  depend  obviously  on  uo> 
f  and  7 . 

Let  y  c  H  be  a  set  of  the  (exact)  solutions  under  consideration.  We 
will  call  the  estimator  g( 7)  a  (M.P.k^  .K^)-proper  if  for  any  uq  e  if 
we  have 
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(2.4) 


k  (M,y)  s  e  s  k  (M,y). 

1  2 

The  estimator  will  be  called  asymptotically  (X.y)-exact  if  for  any  uq  e 
y  and  any  sequence  e  M.,  i  =  1,2...  such  that 

(2.5)  flus(yi)  "  ^01E  0 

as  i  -»  «  we  have 

(2.6)  0(3'i)  ->  1. 

If  0  s  1,  i.e.  a  1,  then  8  will  be  called  iM.y) -upper  estimator 

and  if  0  s  1,  i.e.  k  s  1,  then  8  will  be  called  (jff.y)-lower  estimator. 

2 

We  note  that  k  ,  »c  depend  on  We  cam  also  consider  the  case  that  f 

consists  of  one  element  only,  namely  uq. 

If  Kl>K2  would  be  known,  then  an  upper  ,  respectively,  lower  estimator 

1  1 

can  be  obtained  by  multiplication  of  8  by  —  ,  respectively,  —  . 

l  K2 

An  asymptotically  exact  error  estimator  can  be  easily  constructed  by 
exploiting  the  superconvergence  effect  of  the  derivatives  of  the  finite 
element  solution.  A  systematic  study  of  the  super convergence  phenomena  in  the 
finite  elements  seems  to  have  begun  in  early  1970.  By  now  many  different 
superconvergence  effects  have  been  analyzed.  For  the  reviews  and  references, 
we  refer  to  [20],  [27]  and  [30].  Most  results  on  superconvergence  of  the 
derivatives  for  linear  triangles  are  based  on  [22]  and  [23].  These  ideas  were 
used  very  effectively  in  [19].  See  also  [16]  and  [17].  Superconvergence 
results  hold  only  for  a  very  small  class  ( M,y ),  when  the  meshes  are 
essentially  uniform  and  the  solutions  are  smooth.  In  [1]  it  is  shown  that  for 
quadratic  isoparametric  elements,  the  results  are  of  a  very  similar  nature. 

In  [18]  the  meshes  T^  with  triangles  that  are  obtained  from  triangles  by 
uniformly  subdividing  each  element  of  the  mesh  into  4  congruent 
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subtriangles  are  analyzed.  It  is  shown  that  under  certain  restriction  on  S' 
a  quadratic  interpolation  of  the  solution  obtained  by  linear  elements  yields 
the  same  (asymptotically)  accuracy  in  gradients  as  of  the  quadratic  elements. 
Various  asymptotically  exact  estimates  were  derived  differently  although  the 
superconvergence  is  essentially  present.  See,  e.g.,  [13]  and  references 
there.  The  superconvergence  effects  need  a  smooth  solution  and  regular  mesh. 
This  assumption  can  be  weakened  when  only  superconvergence  on  an  interior 
domain  is  analyzed  [26].  See  also  [21]  for  relevant  results.  The  interior 
estimates  effect  can  be  used  for  proving  asymptotic  exactness  of  the  error 
estimators  on  patchwise  uniform  meshes  as  in  [5]  and  [8]. 

Any  asymptotically  exact  estimator  is  closely  related  to  the  super¬ 
convergence  and  hence  a  very  stringent  assumption  on  the  class  is 

needed.  Nevertheless,  it  is  possible  to  expect  that  the  effectivity  index 
will  be  a  reasonable  one  for  the  class  (jH.^)  which  is  not  too  far  from  the 
one  that  leads  to  the  asymptotic  exactness  of  the  estimator. 

In  contrast  to  the  asymptotic  exactness,  the  properness  holds  for  much 
larger  classes  of  (M,?).  Essentially  it  could  hold  for  any  uq  e  H  defined 
by  (2.1)  and  all  triangulation  which  satisfies  the  minimal  angle  condition. 
Existence  of  ,  k ^  can  be  proven  In  most  cases,  for  example,  using  the 

ideas  in  [2],  [5]  and  [8].  In  Section  3  we  will  give  the  concrete  values  of 

k  ,  k  for  an  estimator. 

1  2 

The  aim  to  establish  the  asymptotic  exactness  of  the  estimator  is  often 
practically  not  important,  because  in  most  applications  the  solution 
is  insufficiently  smooth  and  the  mesh  used  in  practice  is  not  sufficiently 
regular.  On  the  other  hand,  the  mesh  is  often  constructed  adaptively  and  then 
it  is  typical  that  the  error  indicators  are  nearly  equilibrated  ,  i.e.,  that 
we  can  consider  subclass  Mm ,  such  that 
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(2.7) 


M. *  =  jy  e  M  1  (max  rjOf)  /  min  T>(y)l  <  k,  u  e  y\  . 

I  Wey  xey  '  0  ; 

In  addition,  £  could  be  restricted  to  the  functions  which  are  the  solution 
of  the  problem  ?  with  smooth  f,  for  example,  f  e  H1(Q)  and  hence,  uq 
has  special  singularity  in  the  corners  of  the  polygon  £1.  Obviously,  it  is 
practically  important  to  analyze  the  performance  of  an  estimator  for  the  class 
In  Section  5  we  will  address  this  question  in  the  one  dimensional 

setting. 

Although  we  mentioned  only  a  simple  model  problem,  all  the  notions 
obviously  make  sense  in  more  complicated  settings,  say,  for  example,  the 
elasticity  problem. 


3.  ANALYSIS  OF  AN  ESTIMATOR 

For  any  triangular  element  x  e  y,  let  £(x)  be  the  set  of  its  three 
edges.  Further,  let 

V  =  ft  Jt  f 

where  by  |x|  we  denoted  the  area  of  x.  For  each  edge  l  of  the 
triangulation  y  set 

i  ii  3u 
_ s 

a  1  dn 

~  ' 

0  for  (cT 


for  icT 


l 


n. 


where  x  and  x  are  the  two  triangles  sharing  the  edge  l  denoted  by 
u  ‘  H  du  |, 

II  3 II 

x  and  x  ;  n  is  the  normal  to  l  outwards  to  x  .  h? —  denotes  the 

in  out  in  ||  dn  I 

3u 

3 

Jump  of  ^ —  across  the  edge  l;  this  value  is  independent  of  the  choice  of 
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the  direction  of  n 

We  define  the  error  indicator  tj(x): 

(3.D  t)(t)  =  J|x|2(nxf)2  +  |  J  |£|2^J 

le£(z)  J 

where  by  |l|  we  denoted  the  length  of  the  edge  i.  Then  we  have 
Theorem  3.1.  (See  [12]). 

l  h  (T)|f|',T  ]  5  1-0-  "Jt  5 


xe7 


(3.2b) 


[nl/2 

l  »V)|f|^T] 


X6  H 


In  (3.1)  by  |f|  ,  we  denoted  the  H1(x)  seminorm,  p  £  2  integer  and 

j| uq  -  Up||E  is  the  error  of  the  finite  element  solution  using  elements  of 

degree  p.  Constants  and  are  constants  depending  on  the 

minimal  angle  a  of  the  triangulation. 

Remark.  The  expression  (3.2)  could  be  written  in  more  general  form  and 

also  for  general  boundary  conditions,  but  we  will  not  do  it  here  (for 

more,  see  (121).  The  ideas  of  the  proof  of  (3.2b)  are  also  partially 
related  to  [14]  and  [25]. 

Constants  C^  and  C^  can  be  computed  (see  [12])  and  we  get 
(3.3)  0.548(log  p)  sin  ^  s  C£  s  0.813(log  p)  sm  k 


(3.4) 


0.171  sin1/Z  %  s  C'  5  0.290  sin1/2ff 

2  a  2 


In  Table  3.1  we  give  values  of  C'  and  CH  for  some  a  and  p. 

a  a 
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TABLE  3.1.  Values  of  C'  and  CP  . 

a  a 


0 

a 

C' 

a 

cp 

a 

p  =  2 

n 

P  =  6 

p  =  8 

7.5 

0.051 

2.390 

3.306 

3.660 

3.988 

15.0 

0.072 

1.682 

2.341 

2.609 

2.839 

22.0 

0.087 

1.363 

1.918 

2. 156 

2.343 

30.0 

0.099 

1.169 

1.670 

1.895 

2.058 

37.5 

0.  108 

1.035 

1.508 

1.727 

1.876 

45.0 

0.115 

0.939 

1.400 

1.615 

1.757 

52.5 

0. 119 

0.875 

1.334 

1.547 

1.684 

60.0 

0.  121 

0.850 

1.309 

1.522 

1.657 

Let  us  now  analyze  the  bounds  (3.2)  -  (3.4).  The  bound  is  partially 
asymptotic  because  of  the  term  involving  function  f  without  exact 
specification  of  the  constant  C.  Nevertheless,  this  term  is  known  a 
priori  and  is  much  smaller  than  &  with  the  exception  of  the  functions  f, 
which  are  unlike  in  applications.  If  f  is  constant  on  the  triangles,  which 
is  not  uncommon  in  applications,  then  this  term  is  not  present.  More 
important  is  the  term  ||uq  -  u  ||e>  which  expresses  the  influence  of  the 
smootheness  of  uq.  This  term  can  be  made  arbitrarily  small  by  selecting 
sufficiently  large  p.  Nevertheless,  this  influences  the  constant  CP.  From 
(3.3)  and  Table  3.1,  we  see  that  CP  is  increasing  with  p.  Hence,  when  the 
solution  uq  has  the  typical  singularity  in  the  neighborhood  of  the  corner  of 
Q,  the  estimator  will  underestimate  the  true  error,  especially  when  the  mesh 
is  more  or  less  uniform.  If  the  mesh  is  refined  in  the  area  of  the 
singularity  of  uq,  then  lower  values  of  p  can  be  used  and  hence  the 
estimator  will  be  of  higher  quality.  For  more  about  the  properties  of  the 
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p-version  which  are  essential  here,  we  refer  to  [9]  and  references  therein. 

We  note  that  is  independent  of  p,  which  indicates  that  the 

overestimation  of  the  error  by  8  is  not  influenced  by  the  smootheness  of  the 
solution. 

The  dependence  of  the  performance  of  the  error  estimator  on  the  smoothness 

of  the  solution,  as  discussed  above,  is  typical  for  any  estimators  used  in 

practice.  We  remark  that  it  is  possible  to  show  by  the  arguments  in  [2],  [5] 

and  [9]  that  C ^  has  a  bound  independent  of  p.  Nevertheless,  a  realistic 
a 

value  for  this  bound  is  not  available. 

The  expressions  (3.3)  and  (3.4)  indicate  that  the  estimator  deteriorates 

with  the  angle  a.  We  emphasize  that  the  estimator  does  not  take  into 

consideration  any  relation  between  the  shape  of  neighboring  elements. 

Further,  the  estimate  uses  the  "worst"  element  although  possibly  very  few 

elements  of  this  type  are  present  in  the  triangularization. 

Let  us  remark  that  for  the  elasticity  problem  the  constant  is  almost 

proportional  to  sin1/2|^  as  for  Laplace  equation  but  ~  sin 

while  for  the  Laplace  equation  *  sin-1/  .  This  difference  is  strongly 

related  to  the  Korn’s  inequality.  For  more,  see  [12]. 

We  mention  that  the  estimate  (3.2)  is  the  worst  case  which  is  directly 

related  to  the  (-/H.y)  properness.  The  question  arises  whether  the  estimate 

(3.2)  -  (3.4)  is  optimal,  i.e.,  whether  it  is  achievable  for  some  uq  €  ?  and 

J  €  M.  It  is  possible  to  show  (see  [12])  that  for  a  special  choice  of  the 

(uniform)  mesh  of  triangles  with  angles  ^,0,a  =  ^  -  0  and  a  particular 

solution  u  we  get  for  all  0  s  8  s  ^ 
o  2 

(3.5)  0. 24(  tg  0)1/2  6  =  ||uQ  -  us|e. 

Hence,  if  0  <  0  s  ^  ,  0  is  the  minimal  angle  and  we  get  from  (3.5) 
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(3.6) 


KI(g)[s1n|)'/2«  =  !«0-us|l 

where 

0  <  k  s  K  (0)  s  k  ,  0  £  p  *  %  . 

0  1  l  4 

Hence,  with  respect  to  the  angle  a  of  the  triangle,  the  constant  C'  is 

optimal.  For  ^  s  0  <  ^  the  minimal  angles  is  ^  -  0  =  a  and  we  get 

(3.7)  0.24  (ctg  <x)1/2  8  =  |uQ  -  us|E 

or 

yl/2 

K2(a)[Sin?J  ^  =  lUo“USlE 

with  0<k  sK(a)skk,  0  s  a  s  and  hence,  with  respect  to  the  angle 
0  2  2  1  4 

a  also  the  constant  is  an  optimal  one.  Obviously,  it  would  be  desirable 

to  obtain  the  values  of  and  C^,  which  are  sharp  with  respect  to  a 

particular  class  (.M.f). 

Coming  back  to  the  effectivity  index,  we  obviously  have 


We  note  that  the  notion  of  angle  a  is  related  to  the  differential 
equation  under  consideration.  Here  we  have  considered  only  the  case  of 
problem  T  given  in  Section  1. 

In  applications  we  can  use  the  error  estimator  8 based  on  the  error 
indicator  tj  (7)  with 

7)  (t)  =  K  tj(t) 

where  the  constant  K  can  be  selected  arbitrarily.  For  example,  the  constant 
K  can  be  selected  as  the  geometric  average  of  and  or  it  can  be 

chosen  so  that  an  upper  estimate  will  be  obtained  .  Another  possibility  is  to 
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select  K  so  that  the  error  estimator  will  be  asymptotically  exact  if  7  is 
the  mesh  which  satisfies  a  proper  assumption;  in  our  case  it  is  the  mesh  of 
equilateral  triangles  and  K  =  0.26864.  This  value  will  be  used  in  the  next 

section  and  the  estimator  will  be  denoted  by  &* . 

We  note  that  a  correction  factor  K,  which  is  sometimes  used  in  practice, 
cannot  influence  the  dispersion  of  the  effectivity  index  for  class  solutions 
and  meshes. 

4.  SOME  NUMERICAL  EXPERIMENTATION  WITH  THE  ERROR  ESTIMATOR  (3.1) 

In  this  section  we  report  some  computational  results  with  the  estimator  6* 
given  in  the  previous  section.  For  detailed  results  we  refer  to  [11],  The 
aim  of  these  experiments  is  to  see  how  accurately  the  theoretical  results 
presented  in  Section  3  characterize  the  performance  of  the  estimator  in 
various  concrete  cases. 

Example  1..  Consider  Q  =  (0,1)  x  (0,1)  and  the  problem 

-Au  =  f  on  ft 

a  =  0  on  dQ 

and  f  is  chosen  such  that  uq  =  sin  tix  sin  Try.  Further,  consider  the 
regular  3  directional  meshes  that  are  obtained  by  uniform  refinement  of  a 
basic  mesh  with  the  ratio  p  =  i,  n  =  1,2, 4, 8  as  shown  in  Fig.  4.1. 


mesh  ratio  1/1  mesh  ratio  1/2 


Fig.  4.1.  Regular  three  directional  mesh  with  the  ratio  p  =  1,2. 
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TABLE  4.1.  Effectivity  index  of  the  estimator  for  Example  2. 


NUMBER  OF 

ELEMENTS 

MESH  OF  THE  RATIO  p 

1/1 

1.2 

1/4 

1/8 

8 

0.852 

16 

0.993 

32 

1.070 

1.289 

64 

1.183 

1.804 

128 

1.147 

1.503 

256 

1.252 

2.068 

512 

1.180 

1.571 

1024 

1.277 

2.  145 

2048 

1.192 

1.594 

4096 

1 . 287 

2.  169 

UPPER  BOUND 

2.34 

2.87 

3  92 

5.38 

LOWER  BOUND 

0.29 

0.21 

0.  15 

0.  12 

In  Table  4.1  we  report  the  values  of  the  effectivity  index  and  its  bounds 
computed  from  Table  3.1,  for  p  =  2.  From  Table  4.1  we  see  that  as  n  oo, 
the  effectivity  index  converges  to  a  limiting  value  >1.  In  this  particular 
case,  and  asymptotic  expansion  of  0  exists.  The  dependence  on  p  is  not 
the  same  as  for  the  bound, but  relatively  close.  We  see  that  the  performance 
of  the  estimator  deteriorates  with  p  as  expected  from  the  theory.  We  note 
that  0  >  1  in  all  cases,  which, of  course,  does  not  follow  from  the  general 
theory. 

Example  2.  Let  ft  =  (0,1)  x  (0,1)  as  in  Example  1,  and  let  us  consider 
the  problem 


-Au 

=  f 

on  0 

u 

=  0 

for  x  =  0, 1 , 

0 

< 

y  s 

1 

3u 

3y 

=  0 

for  y  =  0,1, 

0 

< 

X  S 

1 

and  f  is  such  that  u 

0 

=  sin  7ix. 

The  meshes 

are 

as 

in 

Example  1,  but 

for  p  s  l  and  p  >  1. 
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This  problem  is  not  exactly  the  one  we  considered  (for  simplicity)  in 
Section  3.  (Nevertheless,  Theorem  3.1  with  the  same  constants  hold  also 


(see  [11]  and  [12] ) . ) 

Table  4.2  shows  the  effectivity  index  of  the  error  estimator. 


TABLE  4.2.  Effectivity  index  of  the  estimator  for  Example  2. 


NUMBERS  OF 

ELEMENTS 

MESH  OF  THE  RATIO  p 

1/1 

1/2 

1/4 

2/1 

4/1 

8 

1.013 

16 

1.564 

0.717 

32 

1.107 

2.258 

0.507 

64 

1.598 

0.783 

128 

1.131 

2.273 

0.553 

256 

1.608 

0.800 

512 

1.137 

2.278 

0.566 

1024 

1.611 

0.804 

2048 

1.139 

2.279 

0.569 

4096 

1.612 

0.806 

UPPER  BOUND 

2.34 

2.87 

3.87 

2.87 

3.87 

LOWER  BOUND 

0.29 

0.21 

0.15 

0.21 

0.  15 

Wee  see  that  the  effectivity  index  converges  as  n  oo,  as  in  Example  1. 

Also,  here  the  asymptotic  expansion  of  0  exists.  The  observed  dependence  on 
p  is  once  more  relatively  close  to  the  theoretical  one.  We  see  that  in  this 
case  the  effectivity  index  can  be  both  larger  and  smaller  than  1  depending 
whether  p  <  1  or  p  >  1. 

Example  3.  Let  us  consider  Q  =  (-1,1)  x  (-1,1)  partitioned 
into  4  square  subdomains  i  =  1,2, 3, 4  as  shown  in  Fig.  4.1 


2 

a(x,y)  =  a 

1 

o(x,y)  =  1 

a 

V 

2. 

ii 

4 

a(x,y)  =  a 

Fig.  4.2.  Coefficient  a  for  the  problem  of  Example  3. 
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Let  us  consider  the  problem 

=  0  on  Q 

=  g  on  3Q 

i  =  1,2, 3, 4,  a(x,y)  =  1  on  u^, 
and  a(x,y)  =  a  on  u  ,  u  .  The  value  a  is  such  that  the  exact  solution 

2  4 

u  in  cj  is 
o 

uQ(x,y)  =  ra(C^  cos(a^)  +  sin(a#>)). 

where  (k,< p)  are  polar  coordinates.  For  a  =  0.5  and  a  =  0.1,  we  get  a  = 
3  +  VS  ~  5.828  and  a  =  166.477...,  respectively.  Although  the  problem 
addressed  in  this  example  is  not  exactly  of  the  type  discussed  in  Sections  1 
and  2,  the  theory  covers  it  too  (see  [11]). 

Let  us  consider  the  uniform  3  directional  meshes  with  the  ratio  p  =  1 
shown  in  Fig.  4.1.  In  Table  4.3  we  give  the  effectivity  index  0. 


TABLE  4.3.  The  effectivity  index  for  Example  3 
and  uniform  mesh. 


NUMBER  OF 
ELEMENTS 

a  =  0.5 

a  =  0.1 

32 

0.658 

0.520 

128 

0.694 

0.512 

512 

0.711 

0.524 

2048 

0.719 

0.536 

9192 

0.722 

We  see  that  when  the  singularity  increases,  the  effectivity  index  is  going 

down  as  predicted  theoretically.  (This  is  the  influence  of  the  term 

II u  -  u  II  in  the  estimate  as  discussed  in  Section  3.  ) 
ii  0  p 11  E 

Let  us  now  use  an  adaptively  constructed  mesh  which  is  refined  in  the 
neighborhood  of  the  origin  where  the  singularity  is  located.  Two  sequences  of 
meshes,  A  and  B,  which  differ  in  the  strengths  of  refinement  were  constructed. 


-7(a(x, y)Vu) 
,  .  3u 

a(x-y)  3K 

where  a(x,y)  is  constant  in  each 
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The  strength  of  the  mesh  refinement  of  the  sequence  B  is  stronger  than  that  of 
A.  In  Fig.  4.3  we  show  two  meshes  of  both  sequences. 


Fig.  4.3.  The  meshes  of  the  sequences  A  and  B. 

Because  of  the  stronger  refinement  in  sequence  B,  a  smaller  value  of  p 
can  be  used  and  hence,  the  effectivity  index  will  be  better  for  sequence  B 


than  for  A. 


In  Table  4.4  we  report  the  effectivlty  index  for  Sequences  A  and  B  of  the 


meshes  for  the  case  a  =  0.5 


TABLE  4.4.  The  effectivity  inces  0  for  Example  3  (a  =  0.5) 
for  the  meshes  with  different  strength  of  refinement. 


SEQUENCE  A 

SEQUENCE  B 

NUMBER  OF 

NUMBER  OF 

ELEMENTS 

0 

ELEMENTS 

0 

32 

0.658 

32 

0.658 

87 

0.712 

64 

0.719 

194 

0.738 

96 

0.787 

370 

0.761 

168 

0.830 

725 

0.787 

244 

0.880 

1297 

0.807 

363 

0.918 

2346 

0.828 

485 

0.960 

604 

0.959 

768 

0.975 

965 

0.995 

We  see  also  that  the  effectivity  index  for  mesh  A  is  better  than  that  for  the 
uniform  mesh  (Table  4.3)  as  predicted  by  the  theory. 

We  have  seen  that  the  theoretical  results  presented  in  Section  3  explain 
very  well  the  computational  data.  For  more  numerical  experiments,  see  [11]. 

Summarizing  the  results  shown  here  and  in  [11],  we  can  make  the  following 
conclusions  about  the  estimator  given  in  Section  3. 

1.  The  estimator  performs  as  expected  from  Theorem  3.1. 

2.  The  effectivity  index  depends  on  the  mesh,  especially  on  the  minimal 
angle  of  the  triangles.  Hence,  avoiding  small  angles  is  preferable.  The 
notion  of  the  angle  depends  on  the  differential  operator. 

3.  The  effectivity  index  is  not  too  sensitive  to  the  topology  of  the  mesh 
uniform  or  not  uniform  (except  for  the  minimal  angle). 

4.  The  effectivity  index  can  be  larger  or  smaller  than  1  depending  on  the 
solution. 
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5.  If  the  solution  has  singular  behavior,  the  quality  of  the  estimator 
deteriorates  and  the  effectivity  index  goes  typically  down.  The  derivation 
cam  be  avoided  by  the  proper  refinement,  for  example,  by  adaptive  procedure, 
as  follows  from  Theorem  3.1. 


5.  ANALYSIS  OF  AN  ESTIMATOR  IN  ONE  DIMENSIONAL  CASE 

As  we  have  said  in  Section  2,  asymptotic  exact  estimators  can  be 
constructed  by  employing  the  superconvergence  phenomena  of  the  recovery  of  the 
gradient  of  the  solution  for  linear  elements  or  higher  derivatives  for 
elements  of  higher  degree.  This  applies  only  for  the  meshes  (essentially 
uniform)  and  solutions  which  produce  the  superconvergence  effects.  Never¬ 
theless,  these  error  estimators  can  be  used  then  for  wider  class  of  meshes  and 
solutions, and  one  can  expect  that  a  K1>K2>  proper  estimator  will  be 
obtained  with  acceptable  values  of  K1>KZ  f°r  a  reasonable  class 
This  is  essentially  the  idea  of  (19]. 

We  will  address  here  a  construction  of  an  error  estimator  in  one  dimension 
which  is  directed  to  meshes  constructed  adaptively.  For  more  details  see 
[10]. 

Let  us  consider  the  model  problem 

(5.1)  -u"  =  f  on  I  =  (0,1) 

u(0)  =  u(l)  =  0 

and  let  f  be  such  that  u  e  H1(I).  Let  7  be  a  partition  of  the  interval 

oo 

I  with  nodes  0  *  x  <  x  <  •••  <  x  =1,  x.  =  (x.  ,,x  J,  h.  =  x.  -  x. 

oi  n  i  i-i  ii  l  l-i 

i  =  l,...,n.  Consider  the  finite  element  method  with  piecewise  linear 
elements  on  J.  In  this  case  we  have 


(5.2) 

Obviously 


Vm'V  =  W- 


u1  ,  ,  is  constant  on  x,  and  hence,  it  is  discontinuous. 
s(x)  i 


Let  us 
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now  construct  a  continuous  piecewise  linear  function  U  so  that  for  i  = 


1 . n-1, 


(5.3) 

where 


u(xi>  =  u;iTl  * 


=  u 


S  1  T 


i+1 


*iKh, 


(5.4) 


'Vi  =  US  I X  "  US  '  T 
1  i+i  i 


is  the  Jump  of  the  derivative  of  ug  at  x^  (obviously  =  1  +  £^),  U(x^) 

depends  only  on  the  values  of  u^  in  the  elements  3^  and  3^+  ,  which 

contains  the  mesh  point  x. .  For  the  definition  of  U(x  )  and  U(x  ),  we  will 

i  on 

extend  the  U(x)  by  an  extension  preserving  the  smoothness  (for  example,  by 
an  antisymmetric  extension). 

We  will  now  define  the  error  indicator 


(5.5) 


tj(t) 


■  u 


2  \W2 
(U(a)-u’s)2  dx 


and  accordingly  the  error  estimator. 


The  problem  Is  now  how  to  select  a. .  Let  A.  =  r- -  and 

l  in 


i+i 


(5.6) 


l.i  1+A, 


(5.7) 


a 

2,1 


1+A 


-l 


The  estimator  6^  based  on  (5.6)  was  proposed  in  [28]  and  the  estimator  8^ 
based  on  (5.7)  was  proposed  in  [5]  and  (91.  Let  us  note  that  both  estimators 
8 i  and  8^  are  the  same  for  uniform  mesh,  i.e.,  for  Ai  =  1.  We  have 


Theorem  5.1  (see  [10]) 

i)  is  an  asymptotically  exact  estimator  for  any  mesh  and  u  e 

H2+C( I ),  e  >  0. 
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li)  &  is  an  asymptotically  exact  estimator  only  for  A  =  1,  i.e.,  for 

the  uniform  mesh.  □ 

Note  that  both  estimators  can  be  understood  as  extension  of  an 
asymptotically  exact  estimator  for  uniform  mesh.  Here  the  estimator  8^ 
asymptotically  exact  for  arbitrary  mesh,  but  similar  result  cannot  be  achieved 
in  two  dimensions. 

Let  us  now  address  the  problem  of  the  optimal  selection  of  a  when  uq 
is  not  sufficiently  smooth  and  Theorem  5.1  is  not  applicable. 

Let  us  assume  that 
i)  u"(x)  >  0  on  I 

ii)  u' (x)  =  u’ I  +  a,  fx-x  ]  ♦  r  (x) 

s'Ti  1 1  i-i  J  > 

2 


where  x 


xi-i+xi 


l-1 
2 


and 


(u-u  )(x)dx 
s 


iii )  (u* 


-u' )|  (x,  )  = 

s  i-i 


aihi 


(1-Li), 


Li  =  o(l) 


ai^i 

(u'-u')|  (x.)  =  — i—  ( 1-R,  ) ,  R.  =  0(1) 

s  'Tj,  i  2  i  l 

2,  3 

r  ah 

I  ((u'-u^)2)dx  =  (1-K2),  =  o(l). 

Ti 

We  remark  that  if  u  is  sufficiently  smooth  L^,  R^,  =  OUO.  The 

assumption  that  u"(x)  does  not  change,  the  sign  in  the  neighborhood  of  the 
singularity  can  be  relaxed. 

We  have 

Theorem  5.2  (see  (10)).  Assume  that  1),  ii),  and  iii)  are  satisfied  and 
that  the  mesh  is  equilibrated  with  respect  to  the  energy  norm  (i.e.,  r^tx)  ~ 
constant),  then 
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(5.8) 


yield  an  asymptotically  exact  estimator. 

Theorem  5.2  shows  that  the  asymptotically  exact  estimator  for  an  unsmooth 


solution  and  adaptively  constructed  meshes  is  different  than  the  estimator 


which  is  asymptotically  exact  for  smooth  solution 


h = 


When  u  is 
o 


smooth  then,  of  course,  A  =  1  +  o(l)  for  the  adaptively  constructed  mesh. 
Let  us  mention  that  the  choice 


(5.9) 


leads  to  the  asymptotically  optimal  estimator  when  the  mesh  is  equilibrated 

with  respect  to  the  norm  II I! l  ’  1  s  s  s  ®,  and  hence,  the  estimator 

s 

is  related  to  the  norm  L  while  the  estimator  6  is  not  a  good  one 

1  2 

for  any  adaptively  constructed  mesh  for  an  unsmooth  solution. 

If  we  would  know  the  values  of  |  (u-u  )  2dx  =  t)  (t),  then  the 

I J  T  S  I  EX 

optimal  error  estimator  will  be  (see  [10])  for 


(5.10) 


i  ^EX  Ti+1 ^  a1/2 
1  +  - i - 7—  A. 

W  1 


w 


T)^) 


we  can  use  these  values  in  (5.10)  and  get  the  values  of  which  will  be 

* 

M2] 

denoted  by  a, 


Let  us  now  present  some  examples. 
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i'* 


'ij 

A 


Example  5.1.  Let  u  =  x*3  -  x  be  the  solution  of  (5.1),  and  let  the  x.  = 

0  /  \  i 

i  =  0,...n.  Further,  let  0.  =  -n — : - r-n -  denote  the  elementary 

1  11“  -Ui»L  (t.) 

2  1 


effectivity  index.  In  Table  5.1  we  report  for  n  =  10  the  value  of  ||u'  - 


s"L  (t, ) 
2  i 


=  error  and  the  elementary  effectivity  indices  0.  for  the  error 


estimator,  based  on  a 


m 


a 
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-[2]  ,  -[2] 
a  and  a 


TABLE  5.1.  The  elementary  effectivity  index. 


T. 

1 

p  =  0.75 

8  =  9.0 

i 

[  1 1 

[2] 

-  [2] 

*  [21 

error 

a 

a 

a 

a 

1 

1 . 99(-3) 

1.290 

1.237 

1.095 

1 .  131 

2 

8.11(-3) 

1.118 

0.881 

1.011 

0.895 

3 

1.42 (-2) 

1.128 

0.861 

0.907 

0.909 

4 

1 . 95 (-2 ) 

1.110 

0.897 

0.929 

0.934 

5 

2. 42 (-2) 

1.084 

0.928 

0.950 

0.957 

6 

3. 86 (-2) 

1.063 

0.948 

0.964 

0.962 

7 

3. 27 (-2) 

1.048 

0.961 

0.973 

0.981 

8 

3. 66 (-2) 

1.037 

0.940 

0.979 

0.971 

9 

4. 04 (-2) 

1.029 

0.976 

0.984 

0.998 

10 

4. 40 (-2) 

0.700 

0.936 

0.986 

0.968 

2.5 


p  =  1.75 


3  = 


i 

error 

a'11 

[2] 

a 

-  [2] 
a 

M2] 

a 

1 

3. 56 (-4) 

3.247 

2.696 

1.036 

1.242 

2 

2. 15 (-3) 

1.656 

1. 373 

1.022 

0.963 

3 

4. 96 (-3) 

1.286 

1.146 

1.005 

1.034 

4 

8.59C-3) 

1.156 

1.076 

1.002 

0.980 

5 

1 . 29 (-2) 

1.097 

1.047 

1.001 

1.019 

6 

1 . 79C-2) 

1.066 

1.031 

1.001 

0.986 

7 

2. 35 (-2) 

1.048 

1.023 

1.001 

1.013 

8 

2. 97 (-2) 

1.036 

1.017 

1.000 

0.989 

9 

3. 64 (-2) 

1.028 

1.013 

1.000 

1.010 

10 

4 . 36 ( -2 ) 

0.  875 

0.  914 

1.000 

0.991 

J 


The  meshes  used  are  approximately  equilibrated  and  hence  the  quality  of 

the  estimator  based  on  a(l1  and  aI2)  are  comparable  and  of  good  quality. 

The  optimally  corrected  estimator  and  the  theoretical  one  (based  on  a[l1  and 
-  [2] 

a  give  a  better  effectivity  index. 

In  Table  5.2  we  show  the  elementary  effectivity  indices  for  the  mesh  which 
is  very  well  equilibrated. 
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TABLE  5.2.  The  elementary  effectivity  index  for  an 
equilibrated  mesh. 


p  = 

=1.75,  p  =  6.5 

(11 

[2] 

error 

a 

a 

1 . 50 ( -2 ) 

1.05414 

0.99399 

1 . 57(-2) 

0.95943 

0.98301 

1 . 57(-2) 

0.99164 

0.99574 

1 . 57(-2) 

0.99608 

0.99809 

1 . 57(-2) 

0.99770 

0.99890 

1 . 57(-2) 

0.99848 

0.99928 

1 . 57 ( -2 ) 

0.99892 

0.99949 

1 . 57(-2) 

0.99920 

0.99962 

1 . 57(-2) 

0.99938 

0.99971 

1 . 57 ( -2 ) 

0.99414 

0.99974 

We  see  from  Table  5.2  that  for  perfectly  equilibrated  mesh,  the  estimator 
[2] 

based  on  a  performs  very  well,  as  expected. 

Let  us  compare  the  performance  of  the  estimator  §  and  &2  based  on 
(5.6  and  (5.7)). 

Example  5.2.  Consider  the  mesh  with  x.  =0.01.  x.  =  0.04.  x.  = 

- K -  j_2  i-1  1 

0.09,  x^+i  =  0.16  and  let  u(x)  =  Re(x^+C^1),  and  let  us  be  interested  in 


0.  for  the  estimators  and  Table  5.3  shows  the  result. 
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TABLE  5.3.  The  elementary  error  indicator  0^  as  function  of  p,q 


=  0.5 


6.  19 
5.21 
2.80 
2.07 
1.69 
1.45 
1.29 
1.  17 
1. 08 
1.02 
0.99 
0.  97 
0.96 
0.97 
0.99 
1.02 
1. 06 
1.11 


p 

0. 

,55 

0. 

65 

0. 

75 

0. 

85 

0. 

95 

1. 

05 

1. 

,  15 

1. 

25 

1. 

35 

1. 

.45 

1. 

55 

1. 

65 

1. 

,75 

1. 

85 

1. 

95 

2. 

.05 

2. 

.  15 

2. 

.25 

1 . 07 
1.03 
1 . 99 
0.97 
0.96 
0.96 
0.96 
0.  97 
0.98 
1.01 
1 . 03 
1.07 
1.  10 
1.  15 
1 .  19 
1.24 
1.30 
1 . 35 


1 . 78 
1 .  68 
1.57 
1.48 
1.39 
1.32 
1 . 26 
1.20 
1 .  15 
1 .  11 
1.07 
1 . 04 
1.02 
1.01 
1.00 
1.00 
1.01 
1 . 02 


4.70 
3.  18 
1.59 
1.  18 
1.01 
0.  94 
0.92 
0.  93 
0.97 


We  see  from  chis  table  that  the  estimator  8  is  better  than  &  for  an 

1  2 

unsmooth  solution,  and  &  is  better  than  6  for  a  smooth  solution,  as 

2  1 

expected,  due  to  the  asymptotic  exactness  of  6 . 

We  can  make  the  following  conclusion: 

1)  The  asymptotically  exact  estimator  for  the  set  of  smooth  solution  does 
not  perform  well  for  an  unsmooth  solution. 

2)  The  asymptotically  exact  estimator  for  the  class  M *  of  equilibrated 
meshes  performs  well  also  for  an  unsmooth  solution.  This  estimator  is  a 
different  one  when  compared  with  the  asymptotically  exact  estimator  for  smooth 
solution. 


3)  For  different  equilibration  criteria  different  estimators  are 
asymptotically  optimal. 

Although  no  two  dimensional  analysis  is  available,  we  expect  that  very 
good  estimators  are  possible  to  construct  for  equilibrated  meshes  as  in  one 
dimension. 
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